###############################################################################-
# Author: Pietryka
# Contact: matthew.pietryka@gmail.com
# Purpose: create Figure 3: Observed vs Permuted (null) Turnout Similarity 
# Notes:
###############################################################################-


# Load Packages  =====================


library(dplyr)
library(tidyr)
library(readr)
library(purrr)
library(ggplot2)
library(glue)

# Load Data ==================================

# load object from "Permutation Tests 2 - Summarize Outcomes.R"
estimates_df <- read_rds("data-files/estimates_df.rds")

# plot preferences.R
source("Plot Preferences.R")



# interpret -----------



estimates_df  %>%
  ungroup()  %>% 
  filter(name == "prop_identical_2016")  %>% 
  mutate(across(where(is.double), round, 2)) %>% 
  t()


estimates_df  %>%
  ungroup()  %>% 
  filter(name == "prop_identical_2018")  %>% 
  mutate(across(where(is.double), round, 2)) %>% 
  t()



estimates_df  %>%
  ungroup()  %>% 
  filter(name == "prop_identical_2020")  %>% 
  mutate(across(where(is.double), round, 2)) %>% 
  t()



estimates_df  %>%
  ungroup()  %>% 
  filter(name == "n_identical")  %>% 
  mutate(across(where(is.double), round, 2)) %>% 
  t()




# Plot functions =====================================




# fixed parameters --------------------
      
the_width <- 2.5
the_height <- 1.5


## plotting function ---------------

plot_density <- function(
    df, # data frame with one row
    outcome_label, # name of similarity measure
    constrained = TRUE, # constrained  (TRUE) or unconstrained (FALSE)
    pct = TRUE # X axis as a percentage (TRUE)?
){
  
  
  if (constrained) {
    sims = df$hist_sims[[1]]
    observed_value = df$observed_estimate[[1]]
    null_mean = df$sim_estimate_mean[[1]]
    p_value = df$p_lab[[1]] 
    fill_color = "grey60"
    diff_estimate = df$diff_estimate[[1]]
    
    
  } else {
    sims = df$hist_sims_unconstrained[[1]]
    observed_value = df$observed_estimate[[1]]
    null_mean = df$sim_estimate_unconstrained_mean[[1]]
    p_value = df$p_lab_unconstrained[[1]] 
    fill_color = "grey85"
    diff_estimate =  df$diff_unconstrained_estimate[[1]] 
  }
 
  diff_lab = ifelse(
    pct, 
    round(diff_estimate * 100, 1) |> paste0("pp"), 
    round(diff_estimate, 2) |> as.character()
    )
  stats_lab = glue("Difference  = {diff_lab}\n{p_value}") 

  
  line_size = 0.8
  line_color = "grey30"
  title_size = 1.4
  
  the_plot <- ggplot(data = sims, aes(x = x, ymin = 0, ymax = y)) +
    format_plot(base_size = 16) +
    xlab(outcome_label) +
    ylab(NULL)  + 
    scale_y_continuous(expand = c(0,  0.06)) +
    scale_x_continuous(labels = scales::number_format(accuracy = 0.1)) +
    geom_ribbon(
      alpha = 0.6,
      fill = fill_color,
      color = line_color,
      linewidth = 0.5
    ) +
    geom_vline(
      xintercept = observed_value,
      color = line_color,
      linewidth = line_size) +
    geom_vline(
      xintercept = null_mean,
      color = line_color,
      linetype = "dashed",
      linewidth = line_size
    ) +
    annotate("text", x = Inf, y = Inf, 
             label = stats_lab, 
             vjust = 1.1, hjust = 1.0,
             size = rel(5)) +
    theme(
      text = element_text(family = the_typeface),
      panel.grid.major = element_blank(),
      panel.grid.minor =  element_blank(),
      plot.title = element_text(
        size = rel(title_size), 
        family = "Arial Narrow Bold"),
      legend.position = "none",
      axis.text.x = element_text(size = rel(1.1)),
      axis.text.y = element_blank()
    )

  if (pct) {
    the_plot + 
      scale_x_continuous(labels = scales::label_percent(accuracy = NULL, scale = 100)) 
  } else {
    the_plot + 
      scale_x_continuous(labels = scales::number_format(accuracy = 0.1)) 
  }

  
}



# Plot Turnout, pooling the genders ==================================




## Outcome: Number Identical -------------------

plot_n_identical_df <- estimates_df  %>% 
  filter(name == "n_identical") 


density_constrained_n_identical <- plot_density(
    df = plot_n_identical_df,
    outcome_label = "Number of Identical Elections,\n2016-2020",
    constrained = TRUE,
    pct = FALSE
    )  +
  expand_limits(x = c(1.3, 2.0))

density_unconstrained_n_identical <- plot_density(
  df = plot_n_identical_df,
  outcome_label = "Number of Identical Elections,\n2016-2020",
  constrained = FALSE,
  pct = FALSE
)  +
  expand_limits(x = c(1.3, 2.0))
  



# interpret
plot_n_identical_df |> 
  select(name, observed_estimate, sim_estimate_mean, diff_estimate, p_lab) |> 
  t()

plot_n_identical_df |> 
  select(name, observed_estimate, sim_estimate_unconstrained_mean, diff_unconstrained_estimate, p_lab_unconstrained  ) |> 
  t()




## Outcome: 2016 -------------------

x_low <- 0.33
x_high <- 0.73
x_lims <- c(x_low, x_high)



plot_2016_df <- estimates_df  %>% 
  filter(name == "prop_identical_2016") 



density_constrained_2016 <- plot_density(
  df = plot_2016_df,
  outcome_label = "Roommates with identical turnout, 2016",
  constrained = TRUE
) +
  expand_limits(x = x_lims)

density_unconstrained_2016 <- plot_density(
  df = plot_2016_df,
  outcome_label = "Roommates with identical turnout, 2016",
  constrained = FALSE
) +
  expand_limits(x = x_lims)





# interpret
plot_2016_df |> 
  select(name, observed_estimate, sim_estimate_mean, diff_estimate, p_lab) |> 
  t()

plot_2016_df |> 
  select(name, observed_estimate, sim_estimate_unconstrained_mean, diff_unconstrained_estimate, p_lab_unconstrained  ) |> 
  t()




## Outcome: 2018 -------------------

plot_2018_df <- estimates_df  %>% 
  filter(name == "prop_identical_2018") 



density_constrained_2018 <- plot_density(
  df = plot_2018_df,
  outcome_label = "Roommates with identical turnout, 2018",
  constrained = TRUE
) +
  expand_limits(x = x_lims)

density_unconstrained_2018 <- plot_density(
  df = plot_2018_df,
  outcome_label = "Roommates with identical turnout, 2018",
  constrained = FALSE
) +
  expand_limits(x = x_lims)




# interpret
plot_2018_df |> 
  select(name, observed_estimate, sim_estimate_mean, diff_estimate, p_lab) |> 
  t()

plot_2018_df |> 
  select(name, observed_estimate, sim_estimate_unconstrained_mean, diff_unconstrained_estimate, p_lab_unconstrained  ) |> 
  t()




## Outcome: 2020 -------------------

plot_2020_df <- estimates_df  %>% 
  filter(name == "prop_identical_2020") 



density_constrained_2020 <- plot_density(
  df = plot_2020_df,
  outcome_label = "Roommates with identical turnout, 2020",
  constrained = TRUE
) +
  expand_limits(x = x_lims)

density_unconstrained_2020 <- plot_density(
  df = plot_2020_df,
  outcome_label = "Roommates with identical turnout, 2020",
  constrained = FALSE
) +
  expand_limits(x = x_lims)



# interpret
plot_2020_df |> 
  select(name, observed_estimate, sim_estimate_mean, diff_estimate, p_lab) |> 
  t()

plot_2020_df |> 
  select(name, observed_estimate, sim_estimate_unconstrained_mean, diff_unconstrained_estimate, p_lab_unconstrained  ) |> 
  t()





# Save ===========================================

windows()
showtext_begin()
showtext_opts(dpi = 600)

save_plot <- function(the_path){
  ragg::agg_png(the_path, 
                width = 180*2, 
                height = 108*2, 
                units = "px"
  )
  
}



save_plot("Results/density_unconstrained_2016.png")
density_unconstrained_2016
dev.off()

save_plot("Results/density_constrained_2016.png")
density_constrained_2016
dev.off()


save_plot("Results/density_unconstrained_2018.png")
density_unconstrained_2018
dev.off()

save_plot("Results/density_constrained_2018.png")
density_constrained_2018
dev.off()

save_plot("Results/density_unconstrained_2020.png")
density_unconstrained_2020
dev.off()

save_plot("Results/density_constrained_2020.png")
density_constrained_2020
dev.off()

save_plot("Results/density_unconstrained_n_identical.png")
density_unconstrained_n_identical
dev.off()


save_plot("Results/density_constrained_n_identical.png")
density_constrained_n_identical
dev.off()

showtext_end()



